package fem2.pre_and_post;

import inf.jlinalg.IVector;

import java.io.IOException;
import java.io.Writer;

import fem2.Element;
import fem2.IntegrationManager;
import fem2.MaterialModel;
import fem2.Mesh;
import fem2.Model;
import fem2.Node;
import fem2.element.StructuralElement;
import fem2.enu.IntegrationType;
import fem2.enu.StrainType;
import fem2.enu.StressType;
import fem2.enu.VariableName;
import fem2.level_set.LevelSet;
import fem2.material.DruckerPragerMM;
import fem2.material.LemaitreDamageMM;
import fem2.material.LemaitreDamageMM2;
import fem2.material.VonMisesMM;

/**
 * class to write post-processing file for GID
 * 
 * @author hbui
 * 
 */
public class GidPostStructural extends GidPost {

	public GidPostStructural(String projectGidDir, String projectName) {
		super(projectGidDir, projectName);
	}

	public GidPostStructural(String projectGidDir, String projectName, String fileName) {
		super(projectGidDir, projectName, fileName);
	}

	/**
	 * write displacement results to nodes. This method is exclusive
	 * 
	 * @param an
	 * @param timeStep
	 * @throws IOException
	 */
	public void writeNodalDisplacementResults(Model m, double timeStep) throws IOException {
		Writer output = createWriter(timeStep);
		output.write("Result ");
		output.write('"' + "DISPLACEMENT" + '"' + " ");
		output.write('"' + "fem2012 analysis" + '"' + " ");
		output.write(timeStep + " ");
		output.write("Vector ");
		output.write("OnNodes\n");
		output.write("Values\n");

		int dim = m.getMesh().getDimension();
		Mesh mesh = m.getMesh();
		for (int i = 0; i < mesh.countNodes(); i++) {
			Node node = mesh.getNode(i);
			double[] u = node.getUHatLocal();
			if (dim == 2) {
				output.write(node.getId() + " " + u[0] + " " + u[1] + " 0\n");
			} else if (dim == 3) {
				output.write(node.getId() + " " + u[0] + " " + u[1] + " " + u[2] + "\n");
			}
		}

		output.write("End Values\n");
		output.close();

		System.out.println("Write nodal displacements completed");
	}

	/**
	 * write nodal level set values to nodes
	 * 
	 * @param m
	 * @param ls
	 * @param timeStep
	 * @throws IOException
	 */
	public void writeNodalLevelSetValues(Model m, LevelSet ls, double timeStep) throws IOException {
		Writer output = createWriter(timeStep);
		output.write("Result ");
		output.write('"' + "LEVEL_SET_" + ls.getName() + '"' + " ");
		output.write('"' + "fem2012 analysis" + '"' + " ");
		output.write(timeStep + " ");
		output.write("Scalar ");
		output.write("OnNodes\n");
		output.write("Values\n");

		Mesh mesh = m.getMesh();
		for (int i = 0; i < mesh.countNodes(); i++) {
			Node node = mesh.getNode(i);
			output.write(node.getId() + " " + ls.valueAt(node) + "\n");
		}

		output.write("End Values\n");
		output.close();

		System.out.println("Write nodal level set " + ls.getName() + " completed");
	}

	/**
	 * write stress at each Gauss point to Gid post file.
	 * 
	 * @param m
	 * @param timeStep
	 * @param checkPrinciple
	 * @throws IOException
	 */
	public void writeGaussPointStressResults(Model m, double timeStep, StressType stressType)
			throws IOException {
		for (String elemType : m.getElementTypes()) {
			/*
			 * extract one element of elemType
			 */
			Element e = null;
			for (int i = 0; i < m.countElements(); i++) {
				if (m.getElement(i).getElemType().equals(elemType)) {
					e = m.getElement(i);
					break;
				}
			}
			/*
			 * extract number of Gauss point
			 */
			if (e == null) {
				throw new Error("element type " + elemType + " does not exist in the model");
			} else {
				int ng;
				if (((StructuralElement) e).getIntegrationType() == IntegrationType.SELECTIVE_REDUCED_INTEGRATION) {
					ng = IntegrationManager.getInstance().getNumIntegrationPoints(
							e.getMeshPart().getType(), e.getIntegrationOrder());
				} else {
					ng = e.getNumIntegrationPoints();
				}
				writeGaussPointStressResults(m, elemType, ng, timeStep, stressType);
			}
		}
		System.out.println("Write " + stressType + " at Gauss point completed");
	}

	/**
	 * write stress at each Gauss point to Gid post file for a specific type of
	 * element and Gauss point
	 * 
	 * @param m
	 * @param elementType
	 * @param ng
	 * @param timeStep
	 * @param checkPrinciple
	 * @throws IOException
	 */
	private void writeGaussPointStressResults(Model m, String elementType, int ng, double timeStep,
			StressType stressType) throws IOException {
		Writer output = createWriter(timeStep);
		int dim = m.getMesh().getDimension();
		String groupName = elementType + "_stress_gp";
		/*
		 * write declaration for GaussPoints block
		 */
		output.write("GaussPoints ");
		output.write('"' + groupName + '"' + " ");
		output.write("ElemType ");
		output.write(elementType + "\n");
		output.write("Number Of Gauss Points: " + ng + "\n");

		output.write("Natural Coordinates: Internal\n");
		// output.write("Natural Coordinates: Given\n");
		// Element dummy = null;
		// for (int i = 0; i < m.countElements(); i++) {
		// StructuralElement e = (StructuralElement) m.getElement(i);
		// if (e.getElemType().equals(elementType)) {
		// dummy = e;
		// break;
		// }
		// }
		// if (dummy == null) {
		// throw new Error("cannot find element of " + elementType);
		// }
		// for (int k = 0; k < ng; k++) {
		// double[] ip = dummy.getIntegrationPoint(k);
		// for (int i = 0; i < ip.length; i++) {
		// output.write(" " + ip[i]);
		// }
		// output.write("\n");
		// }

		output.write("End GaussPoints\n");
		/*
		 * write block of result declaration
		 */
		output.write("Result ");
		output.write('"' + stressType.toString() + "_TENSOR" + '"' + " ");
		output.write('"' + "fem2012 analysis" + '"' + " ");
		output.write(timeStep + " ");
		output.write("Matrix ");
		output.write("OnGaussPoints ");
		output.write('"' + groupName + '"' + "\n");
		output.write("Values\n");
		/*
		 * write data
		 */
		for (int i = 0; i < m.countElements(); i++) {
			StructuralElement e = (StructuralElement) m.getElement(i);
			if (e.getElemType().equals(elementType)) {
				output.write(e.getId() + " ");
				int[] map = getGaussPointMapping(elementType, ng);
				for (int k = 0; k < ng; k++) {
					IVector stress;

					stress = e.getStress(map[k] - 1, stressType);
					// stress = e.getStress(k, stressType);

					/*
					 * write stress
					 */
					if (dim == 2) {
						output.write(stress.get(0) + " " + stress.get(1) + " " + stress.get(2)
								+ "\n");
						// double tmp = map[k];
						// output.write(tmp + " " + tmp + " " + tmp + "\n");
					} else if (dim == 3) {
						output.write(stress.get(0) + " " + stress.get(1) + " " + stress.get(2)
								+ " " + stress.get(3) + " " + stress.get(4) + " " + stress.get(5)
								+ "\n");
					}
				}
			}
		}
		/*
		 * end
		 */
		output.write("End Values\n");

		// /*
		// * start new process: write principle stress
		// */
		// groupName = elementType + "_principle_stress_gp";
		// output.write("GaussPoints ");
		// output.write('"' + groupName + '"' + " ");
		// output.write("ElemType ");
		// output.write(elementType + "\n");
		// output.write("Number Of Gauss Points: " + ng + "\n");
		// output.write("Natural Coordinates: Internal\n");
		// output.write("End GaussPoints\n");
		// /*
		// * write block of result declaration
		// */
		// output.write("Result ");
		// output.write('"' + "CAUCHY_STRESS_TENSOR" + '"' + " ");
		// output.write('"' + "fem2012 analysis" + '"' + " ");
		// output.write(timeStep + " ");
		// // output.write("Vector ");
		// output.write("Matrix ");
		// output.write("OnGaussPoints ");
		// output.write('"' + groupName + '"' + "\n");
		// output.write("Values\n");
		// /*
		// * write data
		// */
		// double o1Max = 0;
		// int eIndex = 0;
		// int gpIndex = 0;
		// for (int i = 0; i < m.countElements(); i++) {
		// StructuralElement e = (StructuralElement) m.getElement(i);
		// if (e.getElemType().equals(elementType)) {
		// output.write((i + 1) + " ");
		// int[] map = getGaussPointMapping(elementType, ng);
		// for (int j = 0; j < ng; j++) {
		// IVector stress = e.getStress(map[j] - 1, StressType.REAL);
		// /*
		// * perform spectral decomposition
		// */
		// double[] pstress = MyLinAlg.spectralDecompose(dim, stress);
		// /*
		// * write principle stress
		// */
		// if (dim == 2) {
		// // output.write(pstress[0] + " " + pstress[1] + " 0\n");
		// output.write(pstress[0] + " " + pstress[1] + " 0\n");
		// } else if (dim == 3) {
		// // output.write(pstress[0] + " " + pstress[1] + " " +
		// // pstress[2] + "\n");
		// output.write(pstress[0] + " " + pstress[1] + " " + pstress[2] +
		// " 0 0 0\n");
		// }
		// if (checkPrinciple) {
		// double o1 = pstress[0];
		// if (o1Max == 0 || o1 > o1Max) {
		// o1Max = o1;
		// eIndex = i;
		// gpIndex = map[j] - 1;
		// }
		// }
		// }
		// }
		// }
		//
		// if (checkPrinciple) {
		// Element emax = m.getElement(eIndex);
		// // StructuralElement semax = (StructuralElement) emax;
		// // IVector sigma = semax.getStress(gpIndex);
		// Debugger.watch("Maximum sigma1 locate at " + emax.getElemType() +
		// " element "
		// + emax.getMeshPart().getId() + " at Gauss point " + gpIndex + " : " +
		// o1Max);
		// // Debugger.watch("sigma =", sigma);
		// // for (int i = 0; i < emax.getNumGaussPoint(); i++) {
		// // Debugger.watch(semax.getStress(i));
		// // }
		// }
		// /*
		// * end
		// */
		// output.write("End Values\n");

		/*
		 * close
		 */
		output.close();
	}

	/**
	 * write strain at each Gauss point to Gid post file.
	 * 
	 * @param m
	 * @param timeStep
	 * @param checkPrinciple
	 * @throws IOException
	 */
	public void writeGaussPointStrainResults(Model m, double timeStep, StrainType strainType)
			throws IOException {
		for (String elemType : m.getElementTypes()) {
			/*
			 * extract one element of elemType
			 */
			Element e = null;
			for (int i = 0; i < m.countElements(); i++) {
				if (m.getElement(i).getElemType().equals(elemType)) {
					e = m.getElement(i);
					break;
				}
			}
			/*
			 * extract number of Gauss point
			 */
			if (e == null) {
				throw new Error("element type " + elemType + " does not exist in the model");
			} else {
				int ng = e.getNumIntegrationPoints();
				if (strainType == StrainType.EQUIVALENT_PLASTIC_STRAIN) {
					writeGaussPointEquivalentPlasticStrainResults(m, elemType, ng, timeStep);
				} else {
					writeGaussPointStrainResults(m, elemType, ng, timeStep, strainType);
				}
			}
		}
		System.out.println("Write " + strainType + " at Gauss point completed");
	}

	/**
	 * write equivalent strain at each Gauss point to Gid post file (refer to
	 * Abdullah)
	 * 
	 * @param m
	 * @param timeStep
	 * @param checkPrinciple
	 * @throws IOException
	 */
	private void writeGaussPointEquivalentPlasticStrainResults(Model m, String elementType, int ng,
			double timeStep) throws IOException {
		Writer output = createWriter(timeStep);
		int dim = m.getMesh().getDimension();
		String groupName = elementType + "_equivalent_strain_gp";
		/*
		 * write declaration for GaussPoints block
		 */
		output.write("GaussPoints ");
		output.write('"' + groupName + '"' + " ");
		output.write("ElemType ");
		output.write(elementType + "\n");
		output.write("Number Of Gauss Points: " + ng + "\n");
		output.write("Natural Coordinates: Internal\n");
		output.write("End GaussPoints\n");
		/*
		 * write block of result declaration
		 */
		output.write("Result ");
		output.write('"' + "EQUIVALENT_PLASTIC_STRAIN" + '"' + " ");
		output.write('"' + "fem2012 analysis" + '"' + " ");
		output.write(timeStep + " ");
		output.write("Scalar ");
		output.write("OnGaussPoints ");
		output.write('"' + groupName + '"' + "\n");
		output.write("Values\n");
		/*
		 * write data
		 */
		for (int i = 0; i < m.countElements(); i++) {
			StructuralElement se = (StructuralElement) m.getElement(i);
			if (se.getElemType().equals(elementType)) {
				output.write(se.getId() + " ");
				int[] map = getGaussPointMapping(elementType, ng);
				for (int k = 0; k < ng; k++) {
					IVector strain;
					int gpMap = map[k] - 1;

					MaterialModel mm = se.getMaterialModel();
					if (mm instanceof VonMisesMM) {
						strain = ((VonMisesMM) mm).getVariables(se).get(gpMap)
								.getVectorVariable(VariableName.PLASTIC_STRAIN);
					} else if (mm instanceof DruckerPragerMM) {
						strain = ((DruckerPragerMM) mm).getVariables(se).get(gpMap)
								.getVectorVariable(VariableName.PLASTIC_STRAIN);
					} else if (mm instanceof LemaitreDamageMM) {
						strain = ((LemaitreDamageMM) mm).getVariables(se).get(gpMap)
								.getVectorVariable(VariableName.PLASTIC_STRAIN);
					} else if (mm instanceof LemaitreDamageMM2) {
						strain = ((LemaitreDamageMM2) mm).getVariables(se).get(gpMap)
								.getVectorVariable(VariableName.PLASTIC_STRAIN);
					}
					/*
					 * TODO: extend for other plastic model here
					 */
					else {
						throw new Error("the equivalent plastic strain does not support for "
								+ mm.getClass().getName() + " material model");
					}

					double eq = 0.0;
					if (dim == 2) {
						eq = strain.get(0) * strain.get(0) + strain.get(1) * strain.get(1)
								+ strain.get(2) * strain.get(2) / 2;
					} else if (dim == 3) {
						eq = strain.get(0)
								* strain.get(0)
								+ strain.get(1)
								* strain.get(1)
								+ strain.get(2)
								* strain.get(2)
								+ (strain.get(3) * strain.get(3) + strain.get(4) * strain.get(4) + strain
										.get(5) * strain.get(5)) / 2;
					}
					eq = Math.sqrt((2.0 / 3) * eq * eq);
					/*
					 * write equivalent strain
					 */
					output.write(eq + "\n");
				}
			}
		}
		/*
		 * end
		 */
		output.write("End Values\n");

		/*
		 * close
		 */
		output.close();
	}

	/**
	 * 
	 * @param m
	 * @param elementType
	 * @param ng
	 * @param timeStep
	 * @throws IOException
	 */
	private void writeGaussPointStrainResults(Model m, String elementType, int ng, double timeStep,
			StrainType strainType) throws IOException {
		Writer output = createWriter(timeStep);
		int dim = m.getMesh().getDimension();
		String groupName = elementType + "_strain_gp";
		/*
		 * write declaration for GaussPoints block
		 */
		output.write("GaussPoints ");
		output.write('"' + groupName + '"' + " ");
		output.write("ElemType ");
		output.write(elementType + "\n");
		output.write("Number Of Gauss Points: " + ng + "\n");
		output.write("Natural Coordinates: Internal\n");
		output.write("End GaussPoints\n");
		/*
		 * write block of result declaration
		 */
		output.write("Result ");
		output.write('"' + strainType.toString() + "_TENSOR" + '"' + " ");
		output.write('"' + "fem2012 analysis" + '"' + " ");
		output.write(timeStep + " ");
		output.write("Matrix ");
		output.write("OnGaussPoints ");
		output.write('"' + groupName + '"' + "\n");
		output.write("Values\n");
		/*
		 * write data
		 */
		for (int i = 0; i < m.countElements(); i++) {
			StructuralElement se = (StructuralElement) m.getElement(i);
			if (se.getElemType().equals(elementType)) {
				output.write(se.getId() + " ");
				int[] map = getGaussPointMapping(elementType, ng);
				for (int k = 0; k < ng; k++) {
					IVector strain;
					int gpMap = map[k] - 1;

					if (strainType == StrainType.SMALL_STRAIN) {
						strain = se.computeStrain(se.getIntegrationPoint(gpMap));
					} else if (strainType == StrainType.PLASTIC_STRAIN) {
						MaterialModel mm = se.getMaterialModel();
						if (mm instanceof VonMisesMM) {
							strain = ((VonMisesMM) mm).getVariables(se).get(gpMap)
									.getVectorVariable(VariableName.PLASTIC_STRAIN);
						} else if (mm instanceof DruckerPragerMM) {
							strain = ((DruckerPragerMM) mm).getVariables(se).get(gpMap)
									.getVectorVariable(VariableName.PLASTIC_STRAIN);
						} else if (mm instanceof LemaitreDamageMM) {
							strain = ((LemaitreDamageMM) mm).getVariables(se).get(gpMap)
									.getVectorVariable(VariableName.PLASTIC_STRAIN);
						} else if (mm instanceof LemaitreDamageMM2) {
							strain = ((LemaitreDamageMM2) mm).getVariables(se).get(gpMap)
									.getVectorVariable(VariableName.PLASTIC_STRAIN);
						}
						/*
						 * TODO: extend for other plastic model here
						 */
						else {
							throw new Error("the " + strainType + " does not support for "
									+ mm.getClass().getName() + " material model");
						}
					}
					/*
					 * TODO: extends for other strain type here
					 */
					else {
						throw new Error("invalid strain type");
					}

					/*
					 * TODO: for 2d model, sigma33 is also available for
					 * PLANE_STRAIN and epsilon33 available for PLANE_STRESS,
					 * even for plastic strain the epsilon33 is always
					 * available. Check this!!!
					 */

					/*
					 * write strain
					 */
					if (dim == 2) {
						output.write(strain.get(0) + " " + strain.get(1) + " " + strain.get(2)
								+ "\n");
					} else if (dim == 3) {
						output.write(strain.get(0) + " " + strain.get(1) + " " + strain.get(2)
								+ " " + strain.get(3) + " " + strain.get(4) + " " + strain.get(5)
								+ "\n");
					}
				}
			}
		}
		/*
		 * end
		 */
		output.write("End Values\n");

		/*
		 * close
		 */
		output.close();
	}
}
